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We consider ray propagation in a waveguide with a designed sound-speed profile perturbed by a range- 
dependent perturbation caused by internal waves in deep ocean environments. The Hamiltonian formalism in 
terms of the action and angle variables is applied to study nonlinear ray dynamics with two sound-channel 
models and three perturbation models: a single-mode perturbation, a random-like sound-speed fluctuations, and 
a mixed perturbation. In the integrable limit without any perturbation, we derive analytical expressions for ray 
arrival times and timefronts at a given range, the main measurable characteristics in field experiments in the 
ocean. In the presence of a single-mode perturbation, ray chaos is shown to arise as a result of overlapping 
nonlinear ray-medium resonances. Poincare maps, plots of variations of the action per a ray cycle length, 
and plots with rays escaping the channel reveal inhomogeneous structure of the underlying phase space with 
remarkable zones of stability where stable coherent ray clusters may be formed. We demonstrate the possibility 
of determining the wavelength of the perturbation mode from the arrival time distribution under conditions of 
ray chaos. It is surprising that coherent ray clusters, consisting of fans of rays which propagate over long ranges 
with close dynamical characteristics, can survive under a random-like multiplicative perturbation modelling 
sound-speed fluctuations caused by a wide spectrum of internal waves. 

PACS numbers: 05.45.Ac; 05.40.Ca; 43.30.+m; 92.10.Vz 



I. INTRODUCTION 

Low-frequency acoustic signals may propagate in the deep ocean to long ranges (up to a few thousands kilometers) due to 
existence of the underwater sound channel which acts as a waveguide confining the sound waves within a restricted water volume 
and preventing their interaction with the lossy ocean bottom jlj. In the ray approximation, the underwater sound propagation can 
be modelled by a Hamiltonian system representing a nonlinear oscillator driven by a weak nonstationary external perturbation. 
A range-independent background sound speed profile plays the role of an unperturbed potential on which a range-dependent 
perturbation of the sound speed along the waveguide, that can be caused by internal waves, mesoscale eddies, ocean fronts or 
something else, is superimposed. 

In the first papers on this topic QQQ], extremal sensitivity of ray trajectories to the initial conditions — ray dynamical chaos 
— has been found in simplified models of the waveguide. In a number of recent publications it has been realized 

that ray chaos should play an important role in interpretating measurements made in the long-range field experiments liolfnl 
which have been designed as a basis for ocean-acoustic tomography JL?! IT3I1 — determining spatio-temporal variations of the 
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hydrological characteristics on the real time scale from acoustical data. The sensitivity of chaotic rays to initial conditions and 
small variations of the environmental parameters causes a smearing of some timefront segments (representing time arrivals in the 
time-depth plane) that has been really observed in the field experiments liollllll . Ray chaos seems to pose restrictions on the ray 
perturbation theory based applications to the tomography. On the other hand, numerical experiments fl^l I15I [hI \y\ fli^l show 
that even at long ranges there exist some stable characteristics of the sound signal which result in remarkably stable segments of 
the timefront — the main measurable characteristic in field experiments used to reconstruct variations in the ocean environment. 
In the recent paper ll^l . maxima of the distribution function of the ray travel time, which lead to clustering of rays, have been 
analytically found with a simplified speed profile corresponding to a quartic oscillator. It has been shown in 1 20] that stable 
fragments of the timefront may correspond to regions of stability in the phase space. Ray stability and instability are strongly 
influenced by the form of the background sound speed profile. 

The ray chaos studies have been especially encouraged by the field experiments [11] where acoustical signals with 75 Hz 
center frequency and 37.5 Hz bandwidth, transmitted near the sound channel axis in the eastern North Pacific Ocean, have been 
recorded with a vertical receiving array between depths of 900 m and 1600 m at a range of 3250 km. The measurements have 
shown a clear contrast between well-resolved earlier portions of the received wavefronts, corresponding to steep rays with large 
values of the action variable, and smearing rear segments of the wavefronts corresponding to near-axial rays with small actions. 

In this paper we study propagation of sound rays in a deep-ocean waveguide with typical sound-speed profiles under internal- 
wave induced single-mode, random-like, and mixed perturbations with the aim to explain and describe peculiarities of ray chaos 
and ray clustering that have been found in natural and numerical experiments. The paper is organized as follows. In SecHTlwe 
give a brief description of the Hamiltonian formalism in terms of the depth-momentum and action-angle canonical variables. In 
Sec. Unlwe design analytically a background sound-speed profile, modelling typical natural deep-ocean profiles. We integrate 
in quadratures the ray equations of motion in the range-independent environment and derive exact expressions for the angle and 
action variables in terms of the depth-momentum variables. Based on the designed profile, we consider two models of sound 



propagation. In Model 1 (Sec. lIII At we exclude from consideration the rays interacting with the ocean surface which cannot 
propagate over large distance because the profile parameters are chosen in such a way that practically all of them interact with 
lossy ocean bottom as well. Shifting the Profile 1 upwards, we obtain Model 2 with rays that may interact with the ocean surface 
without interacting with the bottom. In Sec. llVl we derive analytical expressions for ray arrival times and timefronts at a given 
range with our model range-independent waveguides which should be compared with those in a range-dependent waveguide. 

Section[V]contains results of numerical simulation with Models 1 and 2 in the presence of a single-mode perturbation induced 
by an internal wave. We construct Poincare maps in the polar action-angle variables which show chains of regular islands 
(corresponding to different ray-medium nonlinear resonances) surrounded by a chaotic sea. A new insight into the phase-space 
structure is provided by plots which show by color modulation, respectively, values of variations of the action per ray cycle 
length and values of the range where rays interact with the bottom in terms of initial values of the action and angle variables. 
In the end of this section we demonstrate the possibility of determining the wavelength of the perturbation from arrival time 
distribution under conditions of ray chaos with our model profiles and the Munk canonical one. 

In Sec. lVII we study ray motion under a multiplicative noisy-like perturbation modelling sound-speed fluctuations caused by 
a spectrum of internal waves with flat and decreasing (with the wave number as k~ 2 ) spectral densities. We show that some 
rays may form coherent clusters consisting of fans of rays propagating over long distances with close dynamic characteristics. 
The respective plots of variations of the action are used to clarify a mechanism of appearing coherent clusters in local zones of 
stability in the system's phase space that can survive even under a noisy-like perturbation. The clusterization results in appearing 
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prominent peaks in arrival-time distribution functions and manifests itself in timefronts of arriving signals as sharp strips on a 
smearing background and in plots presenting ray travel time versus starting momentum as "shelf-like segments. 

II. HAMILTONIAN EQUATIONS OF RAY MOTION IN AN UNDERWATER ACOUSTIC WAVEGUIDE 



Consider a two-dimensional underwater acoustic waveguide in the deep ocean with the sound speed c being smooth function 
of depth z and range r. In the geometrical-optics limit, one-way sound ray trajectories satisfy the canonical Hamilton equations 

Q 



dz 
dr 



dp 
~d~r 



with the Hamiltonian 



H 



■\Jn 2 {z, r)- 



(1) 



(2) 



where n(z, r) — cq/c(z, r) is the refractive index, cq is a reference sound speed, p — nsintj) is the analog to mechanical momen- 
tum, and is a ray grazing angle. Only those rays that propagate at comparatively small grazing angles can survive in the ocean 
at long distances, the other ones attenuate rapidly interacting with the lossy ocean bottom. In the paraxial approximation, the 
Hamiltonian can be written in a simple form as a sum of the range-independent and range-dependent parts | 



H = H +H 1 (r) 



with the terms 



Ho 



1 + ^ 



Ac(z) 



<5c(z, r) 



(3) 



(4) 



l Cq Cq 

where Ac(z) — c(z) — cq, 8c(z, r) describes variations of the sound speed along the waveguide. In deriving Eqs. ©, we used 
the condition \n 2 (z, r) — 1| <C 1, that is valid with natural underwater sound channels, and the approximation n 2 (z, r) — 1 ~ 
— 2Ac(z)/co. Moreover, in the paraxial approximation the expression p ~ tan0 is valid. After making the canonical transfor- 
mation from the variables (p, z) to the action-angle variables (/, #), the Hamiltonian may be written in the convenient form 



H = H (I)+H 1 (I,V,r) 



The action variable is defined as the integral [21 



(5) 



(6) 



with z m i n and z max being the depths of the upper and lower ray turning points, respectively. The angle variable is defined as 
follows: 

z 
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where © is the angular frequency of spatial path oscillations, and 

pdz (8) 



is the generating function. 

In a range-independent waveguide the sound-speed profile does not depend on the range r. In such a waveguide the Hamilto- 
nian Ho remains constant along the ray trajectory, and the ray equations in the action-angle variables are trivial 

dl dH d& dH 

Tr = -ld = Q > * = -3r = a ®> (9) 

with the solution 

7 = / , d = i9b + a)(/o)r, (10) 

where 7q = I(r = 0) and #o = &(r — 0) are initial values of the action and angle, respectively. In a range-independent waveguide 
ray trajectories are periodic curves. In a range-dependent waveguide the Hamiltonian equations in terms of the action and angle 
variables take the form 

dl dHi d$ dHi 

d- r = -^ j? =(0+ -dr- 

The action now does not conserve along the ray path. The equations (II 11 are, in general, nonintegrable and are known to have 
chaotic solutions even under a periodic perturbation Hi 0SBG3B 

III. EXACT SOLUTIONS WITH MODEL RANGE-INDEPENDENT WAVEGUIDES 

In this section we study ray nonlinear dynamics in range-independent waveguides with model sound-speed profiles we have 
designed analytically. Our model profiles seem to be attractive by two reasons: they are typical in shape for natural deep- 
ocean background sound-speed profiles and provide analytical solutions to the ray equations including exact expressions for the 
action-angle variables in terms of the depth-momentum variables and analytical ones for timefronts and ray travel times. The 
model profile, hereafter referred as Profile 1 (or Model 1), is depicted in Fig.^ Practically, all the rays, propagating in the 
corresponding waveguide, that interact with the ocean surface interact with the ocean bottom as well. Because of the strong 
attenuation of sound in the bottom, we will exclude such exceptional rays from consideration in numerical simulation. Shifting 
Profile 1 upwards to some distance, as it is shown in Fig.|2] we obtain Profile 2 (or Model 2) with rays that may interact with the 
ocean surface without interacting with the ocean bottom. Both the models will be considered because some characteristics of 
rays, propagating in the respective waveguides, may differ. The Munk canonical profile, widely used in underwater acoustics, is 
shown in Fig. |2J> for comparison. 

A. Model 1 without reflections of rays from the ocean surface 

In Ref. S 

we have introduced a background sound-speed profile, shown in Fig.[0 that models sound propagation through 

the deep ocean 



c(z) = c 



2 



l-^-(i-or-T) 



0<z</z, (12) 
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FIG. 1: Analytic unperturbed sound-speed profile referred as Profile 1 (or Model 1). 



where y = exp(— ah), h is the lower border of the underwater sound channel that is the ocean bottom, a and b are adjusting 
parameters. In simulation with Model 1, we used the following values of the parameters: a = 1.0 km -1 , b = 0.3, h = 2.5 km, 
and cq = 1500 m/s. The depth of the channel axis, where the speed of sound c(z a ) is minimal, is given by 



1 i 2 

Za = - In ■ 



a l + y' 



(13) 



and the parameter b is connected with c(z a ) as follows: 



b = 



1-7 



1- 



The cycle length of the ray path in the channel is given by 



D = 2 



2k 



(14) 



P ayjyb 2 -2E' 



(15) 



where E = Hq + 1 . The Hamilton equations with the range-independent channel d!2t 

dz 



dr 



= P- 



^ = --ab 2 e- az (l + y-2e- az ) 
dr 2 



(16) 



can be solved exactly 



a 2 b 2 



z(r) = -ln- 
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1 +y-2cos(a>r+# ) 



2o) 2 



(17) 
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P( r ) 



(0Qsin((0r+ # ) 



1 + y-fi cos (<*""+#o) 



(18) 



where 0) — 2n/D and #o are the frequency and the initial phase of spatial oscillations of ray path in the channel, respectively. 
We used the short notation in the solutions (1171 and Jl 81 



()(£; = V Ml- 7 ) 2 + f ■ 



The initial phase with a point source, placed at the channel axis, is 



% Q 
j>o = ± — =F arcsin . 

i i + r 



(19) 



(20) 



The unperturbed separatrix is defined by the value E — 0. It is a trajectory that separates propagating rays touching and not 
touching the bottom. Calculating the canonical variables l|6) and with Model 1, we find the action 



7- 



2E 



(21) 



and the angle 



% \ + y-(2y-\E b 2 )e az 

# = ±-=F arcsin — — - — - . 

2 + Q 



The old canonical variables are the following functions of the new ones: 

a 2 b 2 (1 + y-6(7) cos g) 



(22) 



(23) 



^ v ' y a(i+7-eco S1 j)' v 7 



where 



and the frequency of spatial oscillations is given by 

„(/) = 5^-0*7. (26) 

The maximal (at 7 = 0) and minimal (at E = 0) values of the frequency 0)(7) define the minimal and maximal ray cycle lengths, 
respectively 

4-n 2n 

T^min = TT\ \ s s 7) max = : — . (27) 

flfc(l + y) a* a/7 

The derivative dco jdl is known as a parameter characterizing some nonlinear properties of a sound speed profile 

da o 

— = -a 2 . (28) 
dl 

The Hamiltonian can now be written as a function of the action variable only 

, , (\+y)db b 2 (l-y) 2 a 2 , 

Ho(/)= 2 8 -y 7 ■ (29) 
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B. Normal mode amplitudes of the acoustic field in terms of ray quantities 

In this section we derive an exact analytical expression for normal mode amplitudes of the acoustical wave field in the range- 
independent waveguide d!2i in terms of ray variables whose exact solutions we have found above. The connection between 
ray and modal expansions of wave fields in the range-independent environment is well known The normal modes of the 
unperturbed problem satisfy the wave equation 

Ac(z)~ 



1 d 2 Ym , r jz 



Wm = 0, (30) 

CO J 

where k = 2k£1/cq is the wave number in the reference medium with the sound speed cq, Q. is a carrier frequency, and E m = 
1 +Ho(I m ). The eigenfunctions y/ m (z), which represent normal modes in a range-independent waveguide, are supposed to be 
orthogonal and normalized. They constitute a complete set of basic functions in expanding an arbitrary wave field. 

In the Wentzel-Kramers-Brillouin approximation, the eigenvalues of the action variable /,„, corresponding to the m-th mode, 
are determined by the Bohr-Sommerfeld quantization rule 

kl m =m + ^. (31) 

The ra-th eigenfunction \f/,„(z) between its turning points can be represented as follows: 

Vm(z) = V+(z) + y-(z), (32) 



where 



The phase factor is given by 



if£ (z) = A m exp [±i(kG m (z) - k/A)} . (33) 



G m (z) =G(z,/ m ) = J p m {z)dz. (34) 

-'■mill 

The m-th eigenvalue p m with the model profile d!2i can be easily found from Eq. (0} to be 



Pm(z) = ^2E m + b 2 {l~e-^)(e-"^Y)- (35) 

The integral J34b can be calculated exactly 

^ , 7Cl m b{\ + y) , l + y-2e- az 
G m (z) = — H arcsin - 



2 2a Q„ 



2E m . {l + y)b 2 -{2yb 2 -AE m )e- az p,„(z) 
r -^ arCSm ^ (36) 



where Q m — Q(I m ). The amplitude of the m-th mode function is given by the following exact expression: 



^-J *Q+±g± . (37) 
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FIG. 2: Analytic unperturbed sound-speed profile referred as Profile 2 (or Model 2) and the Munk canonical profile. 



C. Model 2 with rays reflecting from the ocean surface 
By shifting Profile 1 (see Eq. H2\ ) upward to a distance d, we get Profile 2 (or Model 2) depicted in Fig.[2]as a solid curve 



c(z) = c 



■y) 



, 0<z<h, 



(38) 



where y = exp [— a(h + d)] , h is the maximal depth of the ocean, cq = c(h), a and b are adjusting parameters. In simulation 
with Model 2 we have used the following values of the parameters: a = 0.5 km -1 , b = 0.6, h — 4.0 km, d — 0.15 km, and 
co = 1535 m/c. In contrary to Model 1, there exist in Model 2 rays which may interact with the ocean surface without interacting 
with the ocean bottom. We will take such rays into consideration because they can propagate to long distances in the ocean. For 
the surface-bounce rays H > H r , where H r is given by 



H r = -\--(\-e- ad ){e- ad -y). 



(39) 



The cycle length of a ray, reflecting from the ocean surface, is the following: 

2 n-A 

D, = ; 

a sjyb 1 - 2E 

where we used the notation 



(40) 



& r = arcsin 

2 



2\ „ud 



l + y-(2y-4E/b 2 )e 



(41) 



In Fig.|5]we show the dependence of the ray cycle length D on the "energy" E = Hq+1. The respective derivative dD/dE has 
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FIG. 3: The cycle length of the ray path D versus the "energy" E for the unperturbed Profile 2. 



a singularity at E r = 1 +H r . As in Model 1, the value E — defines the unperturbed separatrix. We were able to find exact 
expressions for the action 

p(z = 0) b(l + y 1 + y . l + y-2e- ad n-B r 

I = h - — — arcsin 

Ka a \ 4 2% Q % 



2E 



and the angle 



K 



% l + Y-(2Y-4E/b 2 )e< z + d ^ 

$ r - arcsin ■ — 

2 Q 



% l + Y-(2Y-4E/b 2 )e a{z+d ) 

- + arcsin — — — 

2 Q 



, P>0, 
P<0, 



(42) 



(43) 



for the surface-bounce rays. Under reflections, the ray momentum is given by 

p(z = 0) = ^2E + b 2 ( 1 - e- ad ) {e- ad - y) . 
The depth-momentum canonical variables in the range H > H r are the following functions of the action-angle variables: 

a 2 b 2 [l + y+gcosf^^ + tf)) 



(44) 
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a 2 b 2 
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l+ 7 - e COs(^g-t? + 1? r ) 



2co 2 



0<-d<7t. 



(45) 



10 



®esin(^(# + 7T)) 



1+y+gcos U^(# + 7T) 
p(I,*) = { V ' 

l+y-ecos(^-iJ + ^ 
where © is given by the same formula as in Model 1 



-n< ■& <0, 



0< ■& < 7C, 



(46) 



(0 = ayJyb 2 -2E. (47) 

As to normal modes of the unperturbed waveguide 13 8i . they satisfy the respective wave equation J30i with the Bohr- 
Sommerfeld quantization rule 

'. H-_H r . 

(48) 



kl ir , 



W + -, H<H r , 

/r/-i, H > H r . 



At H < H r , the phase factor G m and the amplitude A m (z) of the m-th mode function are given by Eqs. d36l > and 1371 . respectively. 
At H > H r , we get 



„ b / 2E m 



l + y-( 2y -4£ m /^) *M n 
aicsin — h Xr r {E m ) — — 



(i+y)fc 

2a 



l + Y-2e-^ z+ ^ l + Y-2e- 

arcsin arcsin 



ad 



Q„ 



Q„ 



Pm(z = 0) - p m (z) 



(49) 



The amplitude of the w-th mode function at H > H r is 



A m (z) 



l ab{l + Y)-2a 2 I m 

4(71 - & r ) Pm(z) 



(50) 



IV. RAY ARRIVAL TIMES AND TIMEFRONTS IN RANGE-INDEPENDENT AND RANGE-DEPENDENT WAVEGUD)ES 

Internal waves in the ocean induce lateral variations of the sound speed. As a result, the ray cycle length and the ray action are 
not invariants as in range-independent waveguides but vary slowly along the ray path. Even very small variations of the sound 
speed may cause under typical conditions exponential divergence of rays with initially close grazing angles, the phenomenon 
known as ray chaos |Q]. The model of a "frozen" medium is usually adopted, where one may neglect temporal variations in the 
environment and take into account only its spatial variations due to comparatively small propagation time of sound in the ocean. 
Then variations of the speed of sound may be described by the expression 

Sc(z, r) — 5cnns(z)^(z, r), (51) 

where dc^s is the root-mean-square value of sound-speed fluctuations. Following to Refs. we shall describe the 

fluctuations by the simple formula 

5 Crms (z) = ec | e - 2z/B , (52) 

D 
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where e is a measure of the strength of the range-dependent perturbation and B, the termocline depth scale, is chosen to be 1 
km. Throughout the paper, we use the perturbation models with only longitudinal modes of internal waves, i. e., % (z, r) = £, (r). 
Then the perturbed Hamiltonian may be written as follows: 

H = H {I)+eV(I,6)S{r). (53) 

Let us represent the perturbation in the form of the Fourier series over the cyclic variable # 

V(/,0) = i £ V m (I)e im » +c.c. (54) 



= i 



The equations of motion are 



^ = -^£mV m (/) e ^(r)+c.c, (55) 



dr 2 m=\ 



d-d 2% dV Kf s 

^=d +£ ^ W ' (56) 

The function V (I, #) is an analytical one with the Fourier amplitudes exponentially decreasing with increasing the number m. 
With Model 1 and perturbation (I52> . it has the form 

v(i, *) = a-^i° ( ^ ,) \**[i+r-QV)*»*]. (57) 

Methods of the acoustic tomography are actively used for studying spatio-temporal variations in the ocean on the real time 
scale |3 [3 • When the sound waves propagate over long distances, an effective means for monitoring the medium is based on 
the effect of spatial variations of the sound speed on the signal arrival times, one of the main measurable characteristic in long- 
base acoustical experiments. Extensive field measurements, that have been carried out in recent years JloHnl . showed smearing 
of timefront segments in the rear of the sound pulse. Hardly resolvable microfolds in the late-arriving portions of the timefront, 
to be observable in field experiments, can be reasonably explained by the ray's sensitivity to initial conditions. Without internal 
waves the timefront has a smooth folded accordion shape due to refraction as in Fig. |4] In the presence of internal waves, a 
nonuniformity of ray arrivals along the folded fronts appears (see Fig. I16> . Zooming would reveal the presence of microfolds 
along the macroscopic segments of the timefront under consideration. 

In accordance with the Fermat's principle, ray arrival time to a point R along a waveguide is calculated with the help of the 
Lagrangian L 

R R 

t = — I Ldr=— f(p 2 -H)dr. (58) 

CO J Co J 



At sufficiently long ranges, R/D^> 1, the Lagrangian L may be considered as a function of the action 

L(I)=2n^--H (I). (59) 
Following to Eq. (I58> . ray arrival time to the point R along a range-dependent waveguide is given by 

t = -(L(I)}, (60) 
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where (. . .) means an averaging over r. In a range-independent waveguide, arrival times for long-range paths can be simply 
calculated to be 

t^^lfalg^ (61) 

c c V D J 

with the Lagrangian L being an invariant. With a point sound source, all the invariants are functions of the initial value of 
the momentum p(r = 0) = po, i. e., t — t(po). Ray arrival time at the fixed range R is maximal with axial rays (because the 
sound speed is minimal at the channel axis), and decreases in average with increasing po. If t (po) is a monotonic function and 
a waveguide is range-independent, the so-called timefront, which represent ray arrivals in time-depth plane, can be calculated 
explicitly from the equation for the trajectory (see, for example, Eq. d23l for Model 1) with the parameters being functions of 
the Lagrangian L. 

In order to demonstrate it in Model 1 with the range-independent waveguide shown in Fig. [2 we use Eq. i59i to find the action 



I(L) 



the spatial frequency of nonlinear oscillations 



(0(L) = a[b 



i+r 



-aI(L)), 



and the quantity 



Q(L)=2 



(l + Y)aI(L) a 2 I 2 (L) 



(62) 



(63) 



(64) 



b b 2 

as functions of the Lagrangian L. After substituting Eqs. (1601 . J62t — d64t into the ray-trajectory equation (1231 . we get the timefront 
of the sound signal in the waveguide with the sound-speed Profile 1 

2; 2 

a b 



z{t)~ iln 



l + r-Q(t)cos(co{t)R + $o) 



(65) 



a 2co 2 (t) 

The respective plot, presenting ray depths against arrival times at the range R — 1000 km, is shown in Fig.@] We see a typical 
two-folded accordion-like structure due to refraction with positive and negative values of grazing angle. The ray arrivals are 
spread in a smooth and predictable way with the late-arriving portion of the timefront formed by the axial rays. 

In the presence of a perturbation, timefront can be computed approximately as a sum of representing points of sound pulses 



(66) 



where Zi and f, are the depth and arrival time of the 2-th pulse. With the help of Eq. d23l . the depth for each pulse z may be written 
as a function of initial Iq, final //, and mean, (/), values of the action, which may be considered as independent variables under 
conditions of strong chaos 



a 2 b 2 



■In- 



1 + Y- Q(If) cos (©((/) )R + 0o(/ o )) 



(67) 



2(0 2 (I f ) 

Let us consider now distribution of rays over their arrival times, fit , R), at a fixed range R. In a range-independent waveg- 
uide, it is determined by an initial distribution of grazing angles only. The respective distribution function for rays, started at 
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FIG. 4: Timefront for the unperturbed Profile 1: ray depth z versus ray travel time t at the range 1000 km. 



small grazing angles with 0o ~ po (that may propagate over large distances), is fo(po)- Using the condition of the conserved 
normalization 



J f (t,R)dt = J fo(po)dp , 



we get 



/oM)=/o(po) 



dpo 
dt 



With Model 1 we get from Eqs. J6l}, S), O and J2lT> 



fo(t,R) 



Mpo)c Q [b(l+y)-2U(t)] 
RU^b(l+y)U(t)-U 2 (t) ' 



where the short notation is used 



In a range-dependent waveguide, the distribution function of ray arrival times is given by 

dpo(t) 



f(t,R)=AF(t,R)f (po)- 



dt ' 



(68) 



(69) 



(70) 



(71) 



(72) 



where A is a normalization constant. The function F(t, R) describes the effect of the range dependence of the sound speed on 
the distribution of ray arrival times. It is defined, mainly, by the structure of the phase space of the perturbed system 12011 . We 
shall use the function F(t , R) as a convenient tool for analyzing clusterization of rays. 
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It should be noted that the formulas for arrival times (16 11 and timefronts d65l > and J67l > are approximated ones. Comparing 
Fig. |3 plotted using d65i . with numerical simulation shown in Fig.ll2fa. one can see the difference between the two timefronts 
especially in the neighbourhoods of extrema. However, these formulas provide a correct general image of timefronts and give 
a simple analytical connection between measurable ray characterictics and the phase-space ray variables which can be used to 
explain such peculiarities of timefronts as sharp stripes. 



V. RAY CHAOS IN THE PRESENCE OF A SINGLE-MODE PERTURBATION 



A. The phase-space structure 

In this section we consider a single-mode sound-speed perturbation J5 1 1 with 8c tms (z) given by (I52> and the horizontal 
dependence of the internal-wave induced perturbation given by 

4 (r) = cos&r = cos , (73) 
A 

where A is the wavelength of the internal wave. The Hamilton equations take the form 

^- = -L^W lm ^ + C .z., (74) 

l,m 

d -^ = <*+ e ^ y -7T^ + ^ (75) 
dr 2 , dl 

l,m 

Ikr 0o 

where the new phase VP = # 1 is introduced. Ray trajectories are captured in a ray-medium space nonlinear resonance 

m m 

if the condition 

m(0(I) = Ik, (76) 

is satisfied with / and m being integers. This condition can be satisfied at different values of the action variable / re s corresponding 
to resonant tori. Phase oscillations in vicinities of the resonant tori are described by the universal Hamiltonian of nonlinear 

H u =m Q |©;(/ res )| (A/) 2 + £|V/,„|cosm^ , (77) 
where (Qj(I les ) = d(D(I X es)/dI. The width of the resonance in terms of spatial frequency can be approximately estimated as 



A(0=\(o' 1 \AI = 2^e\(0>\V lm , (78) 

where AI is the width of the nonlinear resonance in terms of the action variable. In accordance with Chirikov's criterion l^ . 
global chaos may arise if 

Aco 

TT7 - (79) 

00) 

i. e., if two nonlinear resonances, centered at co and co + 8(0, overlap. Those resonances that overlap slightly form islands in 
the phase space, areas of stable ray motion in a chaotic sea. Nearby the island's borders, one can find the so-called zones of 
stickiness where chaotic trajectories may be localized for long distances r 1 25]. 
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To visualize the structure of the phase space, we construct a Poincare map integrating numerically the ray equations for 
Model 1 with the single-mode perturbation d73l >. Figure[5]demonstrates such a map with the perturbation wavelength X = 10 km 
and a comparatively weak perturbation strength e = 0.0025. It is a two-dimensional slice of the ray motion in a three-dimensional 
space (I x , I y , r mod A ) with I x = I cos # and I y = /sin # being the polar action-angle variables normalized to the separatrix value 
of the action I s (given by Eq. J21> at E = 0) which is a maximal acceptable value of the action. 

A typical (with Hamiltonian systems) picture with stable islands filled with regular trajectories surrounded by a chaotic sea is 
seen in the figure. The chains with 5 and 6 islands correspond to the primary resonances of the first order (I — 1) with m = 5 and 
6, respectively. The chain with 1 1 islands is located between them and corresponds to the second-order resonance with I = 2 and 
m = 11. Even a higher-order resonance with 16 islands (between the (/ = 1 , m = 5) and (/ = 2, m = 1 1) resonances) is seen in 
Fig-13 Concentration of points near island's boundaries indicates sticky trajectories. 



1 




_1 I . : — i-j . 1 

-1 -0.5 0.5 1 




FIG. 5: Poincare map in the normalized polar action-angle variables for Model 1 with the parameters of the periodic perturbation, X = 10 km 
and e = 0.0025. 

Figure|6]shows the Poincare map with rays that may reflect from the ocean surface (Model 2 with e = 0.005 and X = 10 km). 
In difference from Model 1, a stochastic layer appears inside the separatrix loop in a vicinity of the critical value of the action 
/, = I{H r ) (Fig. |6^). We remind that reflection of rays from the ocean surface in Model 2 occurs at H > H r (see Eq. (I39ll). 
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The origin of the localized stochastic layer can be explained as follows. The distance between the resonances of the m-th and 
m + 1-th orders in terms of spatial frequency is equal to 

k " 2 



8(0 



k 
m 



1 



or 

T 



D 



(80) 



and decreases rapidly with increasing D. Since the ray cycle length D has a local maximum at E = E r (Fig. [3}, the resonance 
overlapping in accordance with J79I is maximal near I(E r ), This stochastic layer is isolated from the separatrix by invariant 
curves. As a result, the respective rays are trapped inside the layer forever and their motion is strongly influenced by the fractal 
microstructure of the stochastic layer 11251 \26L 12711 . The fine structure of the phase space is demonstrated in Fig. |6j> where a 
zoom of the region of the stochastic layer near $ = is shown. Chains of microislands corresponding to primary and secondary 
resonances are seen in the figure. It should be noted that an analogous localized stochastic layer has been found with the Munk 
canonical profile Jlfll^ . 
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FIG. 6: (a) Poincare map in the polar normalized action-angle variables for Model 2 with the parameters of the periodic perturbation, A 
10 km and e = 0.005. (b) Zoom of the small region of the stochastic layer indicated in (a). 



Another way to visualize the phase-space structure is provided by the plot that shows by color modulation values of variations 
of the action during the ray cycle length AI in the plane of initial values of the action and angle variables normalized to the 
separatrix value I s and n, respectively. More exactly, AI is a variation of the action between to successive crossings of the line 
= const by a ray. It depends on the initial value of the range ro and may strongly vary in the chaotic regime. It is a distribution 
of variations of the action over the phase space that has important physical meaning. The number of positive variations of the 
action ("hills") and the number of its negative variations ("hollows") are stable characteristics of the system (independent on 
initial conditions) describing the phase oscillations. The respective map for Model 1 with a single-mode perturbation, presented 
in Fig.Q demonstrates an alternating "hills" and "hollows" corresponding to different values of the phase V P. Due to the phase 
dependence, this structure periodically depends on initial values of the range variable r. The "hills" and "hollows" are separated 
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from each other by "zero lines", which correspond to zero variations of the action per cycle length, and may be "stable" and 
"unstable". The "stable" ones transverse the elliptic points of the potential in Eq. d77i . and the "unstable" ones transverse the 
respective hyperbolic points. It should be noted that "hills" and "hollows" may intersect each other. 

Spatial variations of the sound speed along a waveguide may cause the known effect of ray escaping J^l when some rays 
reach the unperturbed separatrix due to diffusion in the action and quit the sound channel. Those rays are supposed to quit the 
channel which interact with the ocean bottom and therefore attenuate rapidly. First of all, steep rays with comparatively small 
arrival times (corresponding to higher modes of the sound field) will escape. The escaping takes place even under an adiabatic 
perturbation |28] but it has peculiarities under the ray chaos conditions. 




e /7i 

FIG. 7: Plot representing variations of the action AI per ray cycle length for Model 1 (h = 10 km and e = 0.005) in the plane of the normalized 
initial values of the action and angle. Bold lines correspond to zero variations of the action per ray cycle length. 
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FIG. 8: Plot representing the values of the range r, where rays interact with the ocean bottom, for the same model and the same parameters as 
in Fig. [7] White color corresponds to those rays that quit the channel during the first ray cycle, i. e. at r < 70 km. 

In the plot, presented for Model 1 in Fig. [8] color modulates the values of the range r where rays interact with the ocean 
bottom, with white color corresponding to rays which quit the waveguide during the first cycle, i. e. at r < 70 km (which is 
a cycle length of a trajectory nearby the separatrix) whereas the black one corresponds to the values r > 1000 km. As in the 
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case with variation of the action, distribution of the values of range for the rays interacting with the ocean bottom is a stable 
characteristic of the system. Topology of the plot is complicated in those areas in the (Io/Is) — (*W^) plane which correspond 
to the stochastic layer. Its patchiness reflects a complicated inhomogeneity of the phase space. We want to stress that channels for 
escaping are formed in the same places where the "hills" are situated in the plot representing variations of the action per ray cycle 
length (Fig.0. The values of the action variable grow rapidly in such channels resulting in increasing trajectory amplitudes. 
One can see in Fig. [8]black spots corresponding to those ranges of the initial conditions of ray trajectories that never quit the 
waveguide. They correspond to the respective islands of regular motion on the Poincare map (Fig.|5j. The patched dark areas 
between the channels with rapidly escaping rays correspond to the initial conditions with long but finite lengths of escaping. The 
angular structure of the distribution of the lengths of escaping in the phase space is an evidence of a non-ergodic ray diffusion. 
In difference from the plot representing variations of the action per ray cycle length, the plot with escaping rays has a manifested 
patched structure due to fractal properties of the phase space typical for open chaotic Hamiltonian systems with weak mixing 
(for a review see I27L . The hierarchy of islands and chains of islands with sticky zones near the island's boundaries, which are 
repeated at all scales of resolution, produces dynamical traps where representing particles may be trapped for a long distance 
(time). It results eventually in anomalous diffusion and power-law distribution functions Iz^l . It is worthwhile to mention that 
escaping of rays is an analog of trapping (escaping) of chaotically advected passive particles in open hydrodynamic flows (see, 
for example, |29]). 



B. Timefront structure under a periodic perturbation 



Let us consider now the role of the ray-medium nonlinear resonance in forming the structure of timefronts of sound signals. 
Due to nonlinear resonance, the ray arrival time of a sound signal along a given ray, captured in a nonlinear resonance near the 
given value of the action / res , tends to its unperturbed value f(L re s) with increasing the distance J^] 

f(Lres) = ^H (81) 

If the width of the resonance is sufficiently large and if there are sufficiently many rays captured in the resonance, the distribution 
function F(t, R) (see Eq. H2\ ) has a pronounced peak near the value t (Lres)- Moreover, it is possible with the help of F(t , R) to 
find the spatial period of a perturbation mode if the arrival time is unambiguously defined by the ray cycle length 1 20] . Inverting 
the resonance condition H6i . we get 

lD ms =mX. (82) 

There can be several peaks of the function F(t,R) with large amplitudes corresponding to resonances with small m at I = 1. 
Thus, if F(t , R) exhibits at least two distinct peaks with comparable amplitudes, we can determine the period of a single-mode 
perturbation as 

A=£>(f res i)-D(r res2 ), (83) 

where f res i and f res 2 are arrival times corresponding to the two peaks. In general, the values of Z)(f res ) can be calculated numeri- 
cally with the help of Eq. d6 1 i . Our model Profiles 1 and 2 admit analytical calculation of the resonant cycle length. In Model 1, 

for example, one can deduce from Eq. ( 16 11 a formula connecting the ray cycle length D with the arrival time t 

-1 

(84) 
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The upper panel in Fig.|9]shows the function F(t, R) corresponding to the timefront computed at R = 1000 km with Profile 
1, the perturbation wavelength X — 10 km, the perturbation strength e = 0.0025 and the other parameters to be specified in the 
preceding section. Results at the range R = 1000 km are represented in this figure by solid lines with the lower axis showing 
the respective values of travel times. Local concentrations of points in the respective timefront (not shown), correspond to two 
sharp peaks of the function F(t, R) at f res i ~ 671.84 s and f res 2 — 672.54 s, with the left one corresponding to the resonance 
(1 = 1, m = 6) and the right one belonging to the resonance (/ = 1, m — 5). Using the formula (18 41 . one can estimate the 
respective resonant cycle lengths, D(t res i) — 60.6 km and D(f res 2) — 50.0 km, and find numerically the perturbation wavelength, 
A ca i ~ 10.6 km. Note that the upper panel in Fig.|9]demonstrates an additional smaller peak at t ~ 672.25 s corresponding to the 
second-order resonance with I = 2 and m = 11 that is manifested on the Poincare map in Fig. [5] as a chain of islands between 
the first-order resonant islands. It has a comparatively small amplitude because the width of this high-order resonance and the 
frequency of phase oscillations are comparatively small. The satellite peak of the primary resonance (I = 1, m = 6) seems to be 
formed by chaotic rays sticking for a long distance to respective resonant islands but quitting this zone somewhere. Note that 
this peak was absent when we have computed the distribution function of ray arrival times at R = 3000 km. 

The upper panel in Fig. ^3 shows the function F(t, R) with Model 1 at the increased value of the perturbation amplitude 
£ = 0.005 for which ray chaos is stronger. Because of a large overlapping of the nonlinear resonances, the peak corresponding to 
the resonance with / = 1 and m~6 has a smaller amplitude than the respective peak in Fig.|9]and disappears at R = 3000 km at 
all. The peak, corresponding to the higher-order resonance with I = 2 and m — 1 1 is absent in Fig.^| The upper panel in Fig. II II 
shows the function F(t, R) with Model 2 at e — 0.005. The left peak corresponds to the resonance (/ = 1, m = 5) which is seen 
on the Poincare map in Fig.|6]whereas the right peak corresponds to the stochastic layer near I(H r ). Since diffusion inside this 
layer is localized the respective rays have close arrival times and form a cluster. 

In order to demonstrate the possibility of determining the wavelength of an internal wave from the ray arrival time distribution 
under conditions of ray chaos not only with our model background profiles, we have computed the timefront and the respective 
function F(t, R) with the Munk canonical background profile (Fig.|2j with a periodic perturbation I23H26 * 

2%r 



c(z, r) = c 



1 + Jl fa - 1+ e^) + <4 e- 2z / B cos — 

B A 



(85) 



where cq = 1500 m/s, jj. = 0.0057, 77 = 2(z — Za)/B is a normalized depth, andB = 1 km. The depth of the channel axis is chosen 
to be Za = 1 km, e = 0.005 and A = 5 km. The results, shown in Fig.^l at the range R = 2000 km, reveal two distinct peaks of the 
function F(t, R), the left one at f res i ~ 1329.9 s corresponds to D(t ms \) ~ 55 km and the right one at f res 2 — 1331.8 s corresponds 
to D(f res 2) — 50 km giving the difference to be equal to the perturbation wavelength A = 5 km. The resonant cycle lengths 
with the Munk profile have been found numerically with the help of Eq. d61> . Generally speaking, determining perturbation 
wavelength is not always possible even in the case of a single-mode perturbation. It is necessary to have at least two chains of 
slightly overlapped primary resonances with 1 = 1, and sufficiently large number of rays should be captured in these resonances. 

Due to stable motion of rays, captured in a nonlinear resonance, their initial and final values of the action are in a comparatively 
narrow interval 81 approximately proportional to the resonance width AI. In accordance with Eq. (1671 . depths of arrivals of 
the resonant rays at a given range R are distributed over narrow intervals. In such a way, double sharp strips (see Fig. I12h). 
corresponding to ray clusters with positive and negative launch angles, appear in the respective timefronts. 
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FIG. 9: Normalized function of distribution of ray arrival times for Model 1 with the periodic perturbation (A = 10 km and e = 0.0025) and 
an imposed multiplicative noise with different values of its strength a. The solid lines represent results at the range 1000 km (see the lower 
axis for travel times), while the dashed lines are computed at the range 3000 km (see the upper axis for travel times). The upper panel shows 
the function with a purely periodic perturbation. 
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FIG. 10: The same as in Fig.|9|but with the parameters of the periodic perturbation, X = 10 km and e = 0.005. 



VI. RAY MOTION IN THE PRESENCE OF A MULTIPLICATIVE NOISE 



A. Ray equations with noise 

Internal waves in the deep ocean are known to have a broadened continuous spectrum of horizontal wavenumbers k which 
may be adequately described by the empirical Garrett-Munk spectrum 1 3 1 ] with the spectral energy density decreasing with 
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FIG. 11: The same as in Fig.[lO|but for Model 2. 

increasing k. The periodic internal-wave induced perturbation of the sound speed to be considered in the preceding section is a 
useful but not realistic approximation. From the theoretical point of view, the sound-speed perturbations due to internal waves 
should be considered as a random function with given statistical characteristics. Under a noisy-like perturbation, there are no 
specified resonances in ray dynamics which, however, is strongly influenced by the presence of a nonlinear sound-speed profile 
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FIG. 12: Timefront and the corresponding function of distribution of ray arrival times for the Munk canonical profile with the parameters of 
the periodic perturbation, A = 5 km and £ = 0.005. 



Let ^ (r) be a stationary stochastic perturbation representing the sound-speed fluctuations caused by internal waves which is 
defined as a spectral decomposition 



;« = / s(k) 



e - ikr dL 



(86) 
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where S(k) = S (&)e^ W and <j)(k) is a random function of the wave number k distributed equally over the interval [0 : 2?r]. The 
perturbation is assumed to be a Gaussian process with normalized first and second moments 

(§(r))=0, {e{r)) = \- (87) 

We assume that the horizontal scale of the internal-wave induced fluctuations is much less than the scale of the range variations 
of the action, and the diffusion approximation can be adopted. The perturbation amplitude depends on the angle variable # and 
is maximal near the upper turning point, # = 0. So, the ray cycle length D gives us a characteristic scale of the range variations 
of the action. In simulation, we realize the process | (r) as a sum of a large number (> 1000) of harmonics distributed in the 
range k € [2^/100 : 2%j\\ km -1 . Two spectral models of the sound-speed fluctuations, So(k) = const and So(k) °< k~ 2 , have 
been used. 

The Hamilton equations of motion 

d? = - £ M^ ^ = C0(I) + £ 1)T^) (88) 

provide the description of sound-ray trajectories through the deep ocean with a broad spectrum of internal waves inducing the 
respective spectrum of the sound-speed fluctuations. Substituting the Fourier decomposition ( I54t in the first equation ( 18 8 1 and 
taking into account that # = co(r)r + i9o, we can write down the variation of the action over the period as 

A/ = e f>V m (/)a,„, (89) 



. . «o D 



a m = —J J S(k)e i{ - m °-V r dkdr +c.c. (90) 

o o 



In order to find (o{r) let us rewrite the second equation J88i as 

fl)(r)-fl)(/o)=e^§(r), (91) 

where Iq is the initial value of the action. Assuming the correlation length of the random process E, (r) to be small as compared 
with the ray cycle length, -C D, we introduce the "fast" variable x — r/r^ which is connected with the angle variable as 
follows: 

x=^. (92) 

The variable co(x) is now treated as a function of x and modelled as a Markovian process with independent increments charac- 
terized by the normal distribution 

< 2 

j7MG ^ V \ 2a 2 

with the variance depending on x and #o 



(94) 



2 J \ dl 
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The main contribution to the integral ( I90> provides the points of a stationary phase given by the wavenumbers k m = mco(x). So 
we get 

Q. m = n{S m {k m {x)]cos[m'd -<l)[k m (x)]]) km ^ = % J S m [k m (x)]cos[nvdo - <j)[k m (x)]]p(k m )dk m , (95) 

where p(k m ) is the respective probability density. Because the perturbation and its derivative dV/dl^) have a sharp 

maximum in a neighbourhood of # = and are approximately zero outside, we can assume p(k m ) ~ v(k m ) at a = a{x = 
— dx)/(Or^). The integral d90l is now given by 

Q m = nSf cos (mt9o - <ftf ) , (96) 
where 0,^ ff and Sf^ are an effective phase and amplitude, respectively. The amplitude is 

= / S (k m )p{k m )dk m , (97) 
o 

where k m = m(0. As a result, we find the variation of the action over the ray period 

M = £7t f) mV m {I)Sf cos (m# - Cf), (98) 

m— 1 

as a sum of "resonant" terms. 

As in the case with a single-mode perturbation, we compute plots which show by color modulation values of variations of 
the action per a ray cycle length, A/, in the plane of the normalized initial values of the action and angle variables. The plot in 
Fig-Elis computed with Model 2 and the fiat spectrum So(k) = const (e = 0.005). To visualize the borders between positive 
("hills") and negative ("hollows") values of A/, the lines with A/ = are bolded in Fig.^] In the range of comparatively small 
values of the action, the first Fourier harmonic V\ in the expansion (I98> is expected to be dominant. Really, only one "hill" is 
present in the figure in the range < Io/L < 0.2. With increasing the action values, the higher-order terms in d98l begin to play 
a more significant role, and the number of "hills" is expected to rise in the respective ranges on the plot representing variations 
of the action per ray cycle length. In Fig. ^] we see two "hills" in the range 0.2 < Iq/I s < 0.4 and three "hills" in the range 
0.4 < /„//, < 0.9. 

Consider now Model 2 with another kind of the perturbation spectrum, So °< k~ 2 , and the same perturbation amplitude e = 
0.005. In difference from the case with the flat spectrum, only one large "hill" now presents in Fig.^^in the whole range of the 
action values. Moreover, the maximal variations of the action are larger as compared with the flat spectrum. One may conclude 
that under the noisy-like perturbation with the spectrum So °< k~ 2 the first harmonic in the series d98l > is a dominant one in all the 
accessible phase space. The dependence A/(#o) is close to a cosine-like one. 

Under a noisy-like perturbation, topology of the plots of variation of the action depends randomly on initial values of the 
time-like variable r, and this dependence is stronger in the case of the flat spectrum, especially in the range of large values of the 
action. It may cause not only phase shifts but even the number of "hills" and "hollows" may vary under varying initial values of 
r. In the case with decreasing spectral density, only smooth shifting of a "hill" along the & axis may occur when varying initial 
values of r. 

B. Coherent ray clusters 

The plots of variations of the action allow us to treat the ray motion as a slow diffusion in the phase space between "hills" and 
"hollows". If these "hills" (or "hollows") are sufficiently large there may arise large fans of rays with close initial conditions 
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FIG. 13: Plot representing variations of the action AI per ray cycle length for Model 2 under the noisy-like perturbation with So(k) = const and 
e = 0.005. Bold lines correspond to zero variations of the action per ray cycle length. 
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FIG. 14: The same as in Fig. ll3l but with So(k) <* k 2 and £ = 0.005. Bold lines correspond to zero variations of the action per ray cycle length. 



preserving close current dynamical characteristics over long distances coherent ray clusters. In the ranges of strong variability 
of the phase space structure, phase correlations decay rapidly resulting in rapidly decaying clusters. So, the length of the phase 
correlations characterizes the stability of a cluster. Similar clusterization may occur in different physical systems (see, for 
example, Therefore, the whole cluster structure may be considered as consisting of statistical and coherent parts. The 

rays, belonging to the statistical part, propagate in the same areas of the phase space with the same value of the Lagrangian (L), 
do not correlate with each other and demonstrate exponential sensitivity to initial conditions. To the contrary, the rays in the 
coherent part do not show sensitive dependence on initial conditions. Two rays with initial values of the momentum po = —0.02 
and po = —0.03 are shown in Fig. ^] in the range interval r € [2000 : 3000] km. The clusterization may influence strongly 
timefronts of sound signals. The prominent stripes visible in the timefront for the stochastic ray simulation (Fig. I16> . which 
belong to ray clusters, resemble the respective strips visible in the timefront fragments for a deterministic perturbation (see 
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FIG. 15: Two ray trajectories with starting momenta py = —0.03 (solid line) and —0.02 (dashed line) for Model 2 under the noisy-like 
perturbation with So(k) <* k~ 2 . 



Fig.ll2h). It should be emphasized that similar strips have been found in the field experiments 




A decoherence and breaking of the respective coherent clusters become prominent with increasing the range (see fuzzy 
segments in Fig. 116b at the range R = 3000 km). It is seen from Fig.ll6B that late-arriving rays are registered not at the channel 
axis (z a ~ 1 km) but rather deeper, at z — 1.5 km. Such a shift of the sound energy down in the depth may be explained as 
follows. The late-arriving signal is formed by a coherent cluster with near-axial rays deflected under propagation from the axial 
value of the action 1 = 0. It follows from Eq. (I67> that rays in this coherent cluster could arrive (at R = 3000 km) at the depth 
different from z a - 

Figure [n] presents the ray travel time t as a function of the starting ray momentum po for Model 2 under the noisy-like 
perturbation with both the spectral models, So(k) = const (a) and So(k) °< k~ 2 (b). All rays are chaotic under a noise perturbation, 
and one might naively expect to see randomly scattered points in the t— po plots. In fact, we see in Fig.^]smooth "shelf-like 
segments alternating with unresolvable structures. Each "shelf" corresponds to a coherent cluster of rays. The "shelves" are 
distributed chaotically over the range of the starting momenta and their positions depend on a specific realization of the random 
process % (r). Comparing between the two spectral models (Figs. I17h and b), we may conclude that the coherent cluster structure 
is more prominent with the spectral model So(k) °< k~ 2 . Such a "shelf-like structure has been found in t(po) plots for a model 
with a single-mode perturbation Jlrjl . with "shelves" to be prescribed to regular islands in the respective phase space. The 
presence of "shelves" may complexify kinetic description of the ray motion with the help of a one-dimensional Fokker-Plank 
equation JT7I I27I1 because the radius of phase correlations is not small in the presence of coherent clusterization. 
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C. Periodic perturbation with a multiplicative noise superimposed 



In the end of this section we consider a perturbation consisting of a periodic dependence of the sound-speed fluctuations on r 
and a multiplicative noise superimposed 

oo 

|(r) = (l-a)cos — + a / S(k)e- ikr dk, 0<a<l, (99) 

where a is the strength of the noisy part. We have computed timefronts of sound signals at different values of a. Fig.|9]shows 
function F(t,R) with Model 1 at the fixed ranges R = 1000 km and R = 3000 km with the perturbation amplitude e = 0.0025 
and spatial period A = 10 km corresponding to the Poincare section in Fig. [5] The solid lines in Figs. I9HTT1 represent results 
at the range 1000 km with the lower axis showing the respective values of travel times, while the dashed lines are computed at 
the range 3000 km with the upper axis for travel times. As it expected, the amplitudes of the prominent peaks, caused by the 
nonlinear resonance with the periodic perturbation, decreases with increasing the values of a. On the other hand, the amplitudes 
of the peaks in the late-arriving signal, caused by the noisy-like perturbation, increases with increasing a. 

Fig-EH s h° ws the function F(t, R) for Model 1 under conditions of more strong chaos at increased value of the perturbation 
amplitude £ = 0.005. All the "deterministic" peaks have comparatively small amplitudes even at a = 0.5. The distribution 
function F(t, R) is shown for Model 2 in Fig. ^2 with the same values of the parameters of perturbation as in the preceding 
figure with Model 1. It is seen that the peaks at t ~ 371.3 s (R = 1000 km) and at t = 2014 s (R = 3000 km), corresponding to 
the primary resonance of the first order (I = 1, m = 5), disappear when the strength of noise reaches the magnitudes of the order 
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FIG. 17: Ray travel time versus starting momentum for Model 2 under noisy-like perturbation with (a) So(k) = const and (b) So °= k 



a >0.5. 

In difference from Model 1, coherent ray clusters in Model 2 appear not only in the late-arriving portion of the signal but, as 
well, in the early arriving portion corresponding to rays reflecting from the ocean surface. The rays with H > H r in Model 2 
are less chaotic than the rays with H < H r . These results show that ray dynamics is strongly influenced by the form of the 
background sound-speed profile. Depending on the form of the background profile, coherent ray clusters may appear in earlier, 
middle and later portions of a timefront. In our opinion, a stability in earlier portions of the wavefront to be measured in the field 
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experiments liolfnl could be explained by peculiarities of the respective background sound-speed profile. 

Under a deterministic perturbation with only a few frequencies, chaoticity of rays is defined mainly by the density of overlap 
of nonlinear resonances characterized by the Chirikov's criterion i79t and is connected with the derivative of the frequency of 
spatial oscillations over the action \cOj\ d28i . However, the Chirikov's criterion is hardly applicable under conditions of a noisy- 
like perturbation with a large number of frequencies [32, In this case we propose to use as a criterion of stochasticity the 
rate of decreasing of Fourier amplitudes in the series ( 1981 . Stochasticity of rays has been shown to become stronger if many 
terms present in the series (I98> . It can be used as a stochasticity criterion for nonlinear systems under a noisy-like perturbation 
in a close analogy with the Chirikov's criterion for deterministic dynamical systems since the rate of decreasing of Fourier 
amplitudes is defined mainly by the dependence of the frequency of spatial oscillations co on the action /. It is a linear function 
with Model 1 ( I26l l. In Model 2 the respective dependence co(I) has a local maximum at /,. ~ I{H r ). Thus, there arise conditions 
for forming coherent ray clusters in Model 2. 



VII. CONCLUSION 



We have treated chaotic and stochastic nonlinear ray dynamics in underwater sound waveguides with longitudinal variations 
of the speed of sound caused by internal oceanic waves. Two models of sound-speed profiles, which are typical in shape for 
deep ocean sound channels, were designed analytically. We were managed to derive with them exact analytical solutions to the 
ray equations of motion without perturbation in terms of the depth-momentum and the action-angle variables and to find exact 
expressions for the frequency of spatial ray oscillations, the timefront of the sound signal at a fixed range and ray travel times. 
Three different kinds of internal-wave induced perturbations have been considered: a single-mode perturbation, a noisy-like 
multiplicative perturbation, and a periodic perturbation with a multiplicative noise superimposed. 

We have found coherent clusters consisting of fans of rays with close dynamical characteristics over long distances and close 
arrival times. It is essential that forming the coherent clusters occurs under different kinds of perturbations, as periodic as noisy- 
like ones. The mechanism of their forming has been found to be connected with existence of specific zones of stability in the 
phase space of the perturbed system under consideration. In the case of a periodic perturbation, these zones appear due to ray- 
medium nonlinear resonances. In the case of a noisy-like multiplicative perturbation, zones of stability appear due to selective 
resonant interactions between different spectral components of the perturbation and harmonics of the unperturbed motion. As 
a result, the phase space has a specific "resonant" topology with local zones of stability. In order to visualize the topology, we 
have used the plots of variations of the action per ray cycle length. 

We proposed a criterion for forming the coherent clusters, namely, the rate of decreasing of the Fourier amplitudes of a 
perturbation written in terms of the canonical action and angle variables. The effect of coherent clusterization depends on the 
horizontal spectrum of the field of internal waves as well. The clusterization becomes more prominent if the spectral density 
decreases rapidly with increasing the wave number k. The clusterization results in forming prominent peaks of functions of 
distribution of arrival times and manifests itself in timefronts of arriving signals as sharp strips on a smearing background formed 
by chaotic rays. It is worthwhile to stress that such strips have been found in timefronts measured in the field experiments Q|. 
The clusterization may cause a redistribution of the acoustic energy over the depth, a stability of early arriving part of the sound 
signal and other effects. It should be taken into account in kinetic modelling of ray dynamics. 

From the standpoint of acoustic tomography of the ocean, the coherent clusterization seems to be a useful property for the 
purpose of determining spatio-temporal variations of the hydrological characteristics on the real time scale under conditions of 
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ray chaos. In a more general context, such a clusterization is interesting from the standpoint of general theory of influence of 
external multiplicative noise on Hamiltonian systems. 
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